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ABSTRACT 

Context. Over 30 planetary systems have been discovered to reside in binary stars. As some of the binary separations are smaller than 
20 astronomical units (AU) the gravitational perturbation of the secondary star has a very strong influence on the planet formation 
process, as it truncates the protoplanetary disk, possibly shortens its lifetime, and stirs up the embedded planetesimals. Due to its small 
semi-major axis (18.5 AU) and relatively large eccentricity e = 0.35 the binary star y Cephei represents a particularly challenging 
example worthy of study in greater detail. 

Aims. In the present study we model the orbital evolution and growth of embedded protoplanetary cores of about 30 earth masses in 
the putative protoplanetary disk surrounding the primary star in the y Cep system. 

Methods. We assume coplanarity of the disk, binary and planet and perform two-dimensional hydrodynamic simulations of embedded 
cores in a protoplanetary disk perturbed by a secondary companion. Before embedding the planet, the equilibrium structure of the 
disk for the observed binary parameters of y Cep is determined. We initiate the embedded planets in the disk on circular orbits with 
different initial distances from the primary. 

Results. The presence of the eccentric secondary star perturbs the disk periodically and generates strong spiral arms at periapse 
which propagate toward the disk centre. The disk perturbations then weaken as the secondary approaches apoapse. The disk also 
becomes slightly eccentric (edisk ~ 0. 1 - 0.15), and displays a slow retrograde precession in the inertial frame. Embedded cores 
interact with the eccentric disk, are periodically disturbed by the strong spiral shocks, and also by the eccentric binary. We find that 
their eccentricity evolution depends primarily on the starting position in the disk. For all initial separations (from 2.5 to 3.5 AU) we 
find inward migration of the cores. For initial semi-major axes ap>2.7, however, we find a strong increase in the planetary eccentricity 
despite the presence of inward migration. Only cores which are initially far from the disk outer edge (ap$2.7AU) have a bounded 
orbital eccentricity which converges, after mass accretion, roughly to the value of the planet observed in the y Cep system. 
Conclusions. Even though a close binary system such as y Cep still presents a challenge to planet formation theory, we have shown 
that under the condition that protoplanetary cores can form at around 2.5 AU, it is possible to evolve and grow such a core to form a 
planet with final configuration similar to that observed. 

Key words, planetary systems: formation - accretion disks hydrodynamics 



1. Introduction 

Among the known 280 extrasolar planets around main sequence 
stars, about 33 are known to reside in binary systems with a wide 
range of orbital separa tions (Eggenbergeret al. 2004; Konackjl 
I2005al: iMugrauer et alj|2005l; iDesidera & Barb ieri 2007). In all 
objects the planet orbits one of the two stars, while the second 
star acts as an external perturber to this syste m. These are the so- 
called S-Type (stellar) orbits (Dvorak' 1986). A truly 'circumbi- 
nary planet' orbiting around a central binary star (a so called 
P-Type orbit) has not yet been detected. 

In most of the systems the binary separation is large enough 
such that the presence of a secondary may not have influenced 
the planet formation process around one of the stars too strongly. 
The subsequent long term evolution of the planet in a binary 
may be shaped by the Kozai effect even for wider binarie s 
(IWu & Murravl llool iTakeda & Raslol 12001 IWuetal.1 llool . 
However, relatively short period binary orbits with binary sep- 
arations of about 20 AU and less have been found in a few sys- 
tems, for example in Gl 86 and notably in y Cephei. This raises 



questions concerning the formation and stability of p lanets in 
such systems. The reality of a planet in HD 188753 ('Konackjl 
[2005ai) . whose hypothetical existence caused some discussions 
about its origin in such a close binary system ( Jang- Condelll 
120071: iPorteeies Zwart & McMiUanll2005llPfahi2005,)r has been 
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put i nto question by newer observations dEggenberger et al.l 
I2007h . 

The mass-period and eccentricity-period distributions of 
extrasolar planets in binaries show interesting trends com- 

pare d with plane ts around single stars, as presented in 

Eggenberger et alj ( 2004.). and more recently in an updated ver- 
sion by Ibesidera & Barbieril ( l2007h . It is argued that for short 
planetary periods ($40 days) the distributions of mass and ec- 
centricities in binaries are significantly different than for longer 
periods. Surprisingly, the closer-in planets in binaries tend to 
have small eccentricities very similar to or even smaller than 
in single stars. Additionally, there seems to be a large fraction 
of high mass planets with small orbital periods in binaries, that 
are not seen in single stars. However, presently the statistics are 
too small to allow definite conclusions. To improve on this defi- 
ciency in the amount of data, several new search programs tai- 
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lored i n particular to p l anets in binary systems h ave been ini- 
tiated ('Konackil l2005bl: iMut erspaug h et"d]|2005h . Considering 
that observationally more and more planets are found to reside in 
binary stars, and may have orbital characteristics different from 
planets orbiting single stars, additional new calculations follow- 
ing the formation and evolution of such systems are urgently 
required. 

Since planets form within protoplanetary disks which are 
tidally perturbed by the presence of a companion star one might 
expect consequences for the planet formation process. The most 
prominent effects are a truncation of the disk and an increase 
in the relative velocities of planetesimals. The influence a sec- 
ondary has on the growth process through collisions from plan- 
etesimals to planeta ry embryos was stud i ed fir st through 3- 
body calculatio n by ^Heppe nheimed (1 19741 1 19781) and later by 
IWhitmire et al.l (11998). Their results implies that for a solar mass 
star the minimum distance to an orbiting companion (perias- 
tron), such that planetesimals at lAU can grow to larger bod- 
ies, must always be larger than about 16 AU. The critical semi- 
major axis scales only weakly with distance from the star In 
the last years these small number studies have been extended 
to full N-body simulations including several hundred planetesi- 
mals. These typicall y utilise a modern N-body integrator adapted 
for binary stars (e.g. [Chambers et aL| l2002). allow for direct col- 
lisions and stickiri^ (iLissauer et al.l 12004). and also gas drag 
(jTurrini et al. 2005|l. For binary parameters such as a Cen ter- 
restrial planet formation from an initial se a of planetesirnals ca n 
proceed similar to the Sun-Jupiter system (ILissauer et al.ll2004l) . 

Another path of study i ng the formation of planets in binary 
stars was taken by Nelson (2000), who studied the interaction of 
two binary stars each of which is surrounded by its own circum- 
stellar disk. The two disks are both modelled hydrodynamically, 
and during each binary o rbit th ey are strongly perturbed at the 
time of periapse. Nelso^ (120001) argues that due to the periodic 
heating of the disk planet formation becomes more unlikely. On 
the one hand, the temperatures in the disk during periapse be- 
come so high that condensation of solid material is reduced, 
while on the other hand the increased temperatures in the spi- 
ral arms hinders t he onset of gravitational instability. However, 
iBossI (119981 12OO6I) has argued that in some cases the presence of 
the secondary may even induce planet formation through gravi- 
tational instability. 

The dynamical evolution of Jovian protoplanets in cir- 
cumbin a ry di sks ( P-type orbits) ha s been investigated by 
iNelsonI (l2003h and iPierens & NelsonI (l2008ah . studies which 
were extended recently to in clude low mass protoplanets 

dPierens & Nelsonll2007ll2008bl) . 

As of today there is only one rather preliminary study (Kley 

I2OOO ) dealing directly with the dynamical evolution of embed- 
ded protoplanets in external binary stellar systems. In that study 
a Jupiter mass planet was evolved in a protoplanetary disk for a 
range of binary separations Ofc,-,, = 50 - lOOAU and a given bi- 
nary eccentricity ei,„ = 0.5. In the binary environment the planet 
experienced a faster growth and migration with respect to the 
unperturbed case. However, the simulations were rather unreal- 
istic in the sense that the planet was not grown in the disk from 
an embryo to its final size, but rather placed 'ready made' into 
an otherwise unperturbed disk. 

In this contribution we would like to improve on these mod- 
els, and in a future publication we will present our results of 
new simulations following the evolution and growth of plan- 
etesimals in the protoplanetary disk in a binary system. In this 
present study we consider the particular system -y Cep because 
this is one of the tighter systems, and thus places rather strong 



Ta ble 1. Orbital data o f the binary y Cep and its planet as given 
by ' Hatzes et"an (l2003l) . which have been used for the simula- 
tions. 
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constraints on the possible planet formation scenarios. For these 
reasons, this system has attracted much attention in past years. 
Several studies looked at the stability and/or t he possibility 
of (additional) habitable planets in the syste m ( Dvorak et al. 
2004; Tu rrini et al.l 1200 4: Haghighi2ou3|2006"; Verr ier & EvansI 
2006). Finally, the possibility to grow protopla nets from very 
small p lanetesimals in y Cep has been studied bv lThebault et al.l 
('20040, who placed km-sized planetesimals in a protoplanetary 
disk including gas drag and gravitational disturbances from the 
secondary and analysed the relative impact velocity. For a suf- 
ficient amount of gas present, the aerodynamic drag between 
gas and particle leads to periastron alignment such that that run- 
away growth is possible inside of 2.5 AU from the primary star. 
A more extended study of this process for various dyna mical 
parameter of the binary has been perfo rmed by Thebau lt et al 



parameter or the bmary rias been perio rme d by Itiebau lt et al 
(2006). Recently , Kiev & NelsonI (l2007l) and lPaardekooper et al 



(2008) improved again on these models by considering the evo- 
lution of a swarm of planetesimals embedded in a hydrodynam- 
ically evolving disk. 

In the pr esent work we focus o n y Cep and we assume, along 
the work of iThebault et alj (|2004|) . that protoplanetary cores of 
a few earth masses have successfully formed inside the disk and 
follow the subsequent dynamical evolution and mass growth. 
We focus here on the later phase of gas accretion only and 
will present additional new res ults on planetesimal formation 
in a forthcoming paper (see also Kiev & Nelsonll2007 ^. We per- 
form full numerical, hydrodynamical evolutions of protoplane- 
tary disks for binary parameters typical for y Cep, analyse the 
equilibrium configuration of the disk, and study the subsequent 
evolution of embedded protoplanets. 

This paper is organised as follows: in Sect.|2]we describe the 
model set-up, and in Sect. |3] we present results of disk evolu- 
tion due to the binary but in the absence of the embedded planet. 
In Sect. |4]we present results of simulations with planets embed- 
ded in a protoplanetary disk perturbed by a binary companion. 
Finally, we draw our conclusions in Sect. |5] 



2. Model Setup 

In modelling the putative early evolution of the y Cep system, we 
assume that the protoplanetary disk is in an S-type configuration 
and orbits around the more massive primary star. To clarify the 
primary effects of planetary growth and evolution in binaries, 
we avoid complex three-dimensional effects and assume that the 
whole system is coplanar, i.e. the orbits of the disk, the binary 
and the planet are all lying in one plane and the angular mo- 
menta are all aligned. To further simplify the calculations we 
treat the disk in a two-dimensional approximation, neglect any 
vertical extension and work with the vertically integrated equa- 
tions. This approximation is often made for accretion disks and 
applies whenever the vertical thickness H (pressure scale height) 
of the disk is small compared to the radial distance r from the 
central star, i.e. h = H/r « 1. 

We use cylindrical coordinates (r, (f, z), where the the origin 
of the coordinate system lies in the centre of the primary star at 
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r - 0,z - 0. The evolution of all objects (stars, disk and planet) 
takes place in the z = plane. As the stars orbit each other, the 
coordinate system is non-inertial and the indirect terms have to 
be included in the equations of motion. 

Accretion disks are most-likely driven by magneto- 
hydrodynamical turbulence which we model here in a simplified 
manner through a viscous stress tensor model. For the kinematic 
viscosity coefficient v we assume a standard a-prescription with 
V - aCsH, where a is a constant, and Cs is the local sound speed 
given by Csir) = h Ukei,{r) where ui,ep{r) is the Keplerian velocity 
of the accretion disk around the primary star. 

The basic hydrodynamic equations (mass and momentum 
conservation) describing the time evolution of such a viscous 
two-dimensional disk with embedded p lanets have been stated 
frequently and are not repeated here (see' Klevlll999l) . 

The dynamical effects of the companion star will lead to a 
disk structure which is different from that which exists around 
a single star. To account for this restructuring of the disk we 
first model the disk evolution in the binary system alone without 
any embedded planet. As will be shown below this initialisa- 
tion phase takes about 100 binary orbits. Only after this time the 
disk has relaxed (on average) to its new equilibrium and then the 
planets are inserted into the models. 

2. 1 . Initial Conditions 

The two-dimensional (r - ip) computational domain consists of 
a complete ring {(pmin = 0, (pmax - 27r) of the protoplanetary disk 
that radially extends from rmin - 0.50AU and rmax - 8.0AU. The 
computational domain is covered by 300 x 300 {Nr x A^^) grid 
cells, that are spaced equidistant in both radius and azimuth. 

The initial hydrodynamic structure of the disk (density, tem- 
perature, velocity) is axisymmetric with respect to the loca- 
tion of the primary star The surface density has the profile 
2(r) = Sor"'''^ over the entire radial domain, where Eo is 
chosen such that the total mass in the computational domain 
(within r„„„ and r„a^-) equals 1.75 x IO^-'Mq which implies 
So = 1.89 X IO^^Mq/AU^. The initial velocity is pure Keplerian 
rotation (ur - 0, = Q.xf - (GM^/ry^-), and the tempera- 
ture stratification is always given by T{r) cc r"' which follows 
from the assumed constant vertical height h - Hjr - 0.05. For 
these locally isothermal models the temperature profile is left 
unchanged at its initial state throughout the computations, i.e. 
no energy equation is solved. 



2.2. Boundary conditions 

The boundary conditions are chosen such that material may es- 
cape through the radial boundaries. At the outer boundary (r,„a^) 
we impose a so called zero-gradient outflow condition, where the 
velocity is extrapolated when directed outward and set to zero 
when directed inward in the last active gridcell inside r^ax- For 
all other variables zero gradients are assumed. During periastron 
when large spiral arms may extend beyond r,„ax this condition 
will allow material to leave the system freely and will help to 
reduce numerical artifacts. At the inner boundary we set a vis- 
cous outflow condition where the material may flow through r„„„ 
with the local (azimuthally averaged) viscous inflow referring to 
an accretion disk in equilibrium, 

3 . 

^r(^mm) — ~ 7 (1) 
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used at r„,a j but found that in such a case the induced elliptic- 
ity of the disk (by the eccentric binary) will cause strong outflow 
of disk material in its periastron phase which in turn produces a 
large eccentric inner hole of the disk which is a numerical arte- 
fact . Our chosen inner boundary allows for accretion onto the 
star and additionally mimics the expected nearly circular flow in 
the inner regions of the disk close to the primary. No matter is 
allowed to flow back into the system and the mass of the disk 
will slowly decUne. 

2.3. Physical parameters 

In the majority of our our numerical simulations we choose to 
model a specific system where the orbital elements of the binary 
have been chosen to match the system y Cep quite closely. The 
data for this system have been taken from Hatzes et al. (2003!) 
which do not include the mo st recent observational improve- 
ments (iNeuhauser et al.l [20071) . These newest refinements pri- 
marily concern the mass of the primary and should not alter 
our main conclusions at all. The present studies are supposed 
to demonstrate the principle physical effects occurring in close 
binaries rather than trying to achieve a perfect match with all 
the observations of this system. For the runs we choose a bi- 
nary with Ml = 1.59Mo, M2 = O.38M0, ati„ = 18.5 AU and 
Sbin - 0.36 (see Table [T]!, which translates into a binary period 
of P = 56.7 yrs. For the viscosity we use a typical value for the 
effective viscosity in a protoplanetary disk, i.e. o- = 5 x 10"^, and 
a locally isothermal disk with h - H/r = 0.05. 

As mentioned above, the initial total disk mass in the system 
is 1.75 X 10^^ Mq. Due to the presence of the secondary the disk 
will be truncated to a smaller outer radius and the radial density 
stratification will adjust to a new equilibrium configuration. At 
the same time mass will be lost from the system due to possible 
outflow through the (open) boundaries. To ensure a somewhat 
uniform initial configuration for the planets we rescale the disk 
mass to a given mass before inserting the planets. 

2.4. A few remarks on numerical issues 

To substantiate o ur results we use two different cod es for our cal- 
culations, RH2D ('KleV'1999', 1989) and NIRVANA (" Nelson et all 
2000; Ziegler & Yorke 1997), both of which have been utiUsed 
frequently in embedded planet problems. The numerical method 
used in both codes is a staggered mesh, spatially second order 
finite difference method based, where advectio n is based on the 
mono tonic transport algorithm ("v an Leer! 19771) . Due to operator- 
splitting the codes are semi-second order in time. The compu- 
tational details of RH2D which can be used in different co ordi- 
nate systems have been described in genera l in KlevI (1 19891) . and 
specifically for planet calculations in lKlevl(ll999l) . The algorith- 
mic details of the NIRVANA code have been described in Zieglga 
(1998). In all of our simulations we utilise a non-rotating coor- 
dinate frame centred on the primary. 

In calculating the gravitational effect of the planet we use a 
smoothed potential of the form 



Of 
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where s is the distance from the planet to the (grid-)point un- 
der consideration. For the smoothing length of the potential we 
choose £ - 0.6/?Hiih where 



where v denotes the azimuthally averaged viscosity at r,„in- We 
experimented with an outflow inner boundary similar to that 



^Hiii = rp 



(3) 
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Here rp denotes the momentary distance of the planet from the 
star, i.e. the size of the Hill radius scales with the distance from 
the star and not the semi-major axis. This distinction is important 
for calculating the torques in the case of eccentric orbits. 

The smoothing length e is proportional to the local scale 
height H of the disk to account for the vertical extent of real 
disks. The chosen value of 0.6/?Hiii gives a reasonably good 
agreement between the linear torques of a three-dimensional 
disk and the two-di mensional flat-disk approximation used here 
dTanaka et al.l2002l) . In a recent study of the dynamical evolution 
of eccentric low mass planets we have found a good agreement 
for the migration and e ccentricity da mping rate in 2D and 3D 
disks for this value of e dCresswell et al. 2007). 

To model possible mass accretion onto the protoplanet we 
take out mass within a given radius race around the planet at a 
prescribed rate. Whenever the centre of a gridcell is closer than 
face = 1 I'^Rhui we reduce the density in that cell by a factor 1 - 
face Af, where At the the actual time-step of the computation and 
face is a model dependent reduction factor. The rate is doubled 
whenever the distance is smaller than 1/4 rnui- For details see 
[Kley (1999). 

The viscous terms, including all necessary tensor compo- 
nents, are treated explicitly. To ensure stability in the outer parts 
of the disk where we expect stronger gradients in the flow due 
to the perturbations by the secondary, an artificial bulk viscosity 
has been added, with a coefficient Can - 1 0. For the detailed 
formulation of the viscosity related issues and tests, see 
([T999.) . 

To avoid excessively low densities in the outer part of the 
computational domain we have found it preferable to work with 
a density floor, where the density cannot fall below a specified 
minimum value Zmin- For our purpose we use a value of Smin - 
10"^ in dimensionless units, where the initial density is of (9(1) 
atr = 1.0. 

When calculating the gravitational forces of the disk acting 
on the planet numerically for a given grid resolution a cutofi" 
radius r, is necessary in addition to the smoothing of the poten- 
tial. This cutoff' reduces the contribution of material in the direct 
vicinity of the planet which lies inside the planetary Roche-lobe. 
This material is bound to the planet and typically cannot exert 
any forces on it. Here we use a cutoff of r, - Q.^Rhui for most 
runs which is smoothed at r, using a Fermi-function with a width 
of O.lr,. The role of the cutoff radius r, in estimatin g the torque 
on the planet has been analysed in more detail bv ICrida et al.l 
(12008). For the range Q.^Rnm < r, < l-ORnm the obtained re- 
sults agree very well with each other, see also Sect. 14. 1 I below. 

The planet evolution and the orbital motion of the binary pro- 
ceed on time scales much longer than the typical dynamical time 
scale of the disk, and consequently very many time-steps need 
to be calculated just in one simulation. To allow for parameter 
studies we have increased the performance of the runs by im- 
plementing the FARGO-algorithm in our codes, which is partic- 
ularly designed to model differentially rotating flows efficiently 
(Masset 2000). For our chosen radial range and grid resolution 
we find a speed-up factor of about 7.5 over the standard case. 
Then, applying a Courant number of 0.75, still about 160,000 
time-steps are required for only 10 binary orbits for our radial 
range and a standard resolution of 300 x 300. We checked the ac- 
curacy of the FARGO implementation against several test cases. 

3. Disk evolution without an embedded planet 

When studying the formation of planets in a protoplanetary disk 
in the presence of a secondary star it is necessary and illustrative 
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Fig. 1. Time evolution of the total disk mass in units of the Mq 
vs. time in yrs for no embedded planet. Note that one binary orbit 
refers to about 57 yrs. 



to first investigate the structure and evolution of the accretion 
disk subject to the influence of the binary system only. Hence, 
in this section we follow the evolution of a disk around the pri- 
mary which is perturbed by the secondary without an embedded 
planet. As the disk evolves slowly from its initial axisymmetric 
state into a new quasi-equilibrium which is quite different from 
the unperturbed initial state, we have first to relax the system into 
this new equilibrium, before adding a planetary embryo at a later 
time. In the following we describe in more detail this restructur- 
ing process of the disk and the change in the binary elements due 
the finite disk mass. 



Fig. 2. Grayscale plot of the two-dimensional density distribu- 
tion of the circumstellar disk around the primary at two different 
orbital phases of the binary. Left shortly after apocentre at about 
20 binary orbits, and Right shortly after closest approach (peri- 
centre). 



3. 1 . The structure of the disk 

The presence of an eccentric secondary star leads to a strong pe- 
riodic disturbance of the disk whenever it is at periastron. Then 
two strong tidal spiral arms (shock waves) are excited in the disk 
which wind themselves all the way down to the centre. At the 
outer edge of the disk the perturbation by the companion star 
carries material beyond the outer boundary of the computational 
domain which is then lost from the system. During the very 
first periapse phase about 30% of the initial disk mass is lost 
from the system, and very little during the following passages, 
as show in Fig.[T] In between the periapses the disk settles down 
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Fig. 3. Radial dependence of the azimuthally averaged surface 
density (Top) and the eccentricity (Bottom) of the circumstel- 
lar disk around the primary in the presence of the secondary at 
different times. Time is given units of the binary orbit, radial 
distance in AU, and the density in dimensionless units. 



and becomes more circular again. This effect is illustrated in the 
Fig. |2]where we display the surface density S of the disk in gray 
scale at 2 different times in the early evolution of the disk. In 
the left panel the nearly circular state at the binary apoastron 
is displayed while the right panel shows the situation just after 
periastron. Th ese findings are very similarly to those found by 
iNelsonI ( |2000|) in his study of the binary system L155 1 . 

Already the very first close approaches with the binary lead 
to a clear truncation of the disk at around 5 AU, as is visible in 
left panel of Fig.[3]for the curve at f = 10 binary orbits. Slowly 
the whole disk structure rearranges and equilibrates at around 
t = 50 where it shows a much steeper (nearly linear) density 
slope than in the initial state. The time scale for this equilibra- 
tion process depends on the magnitude of the disk viscosity and 
shock dissipation. For relatively cool disks and small viscosities 
we expect that the equilibrium profile is primarily determined 



by the same processes, i.e. gravitational action of the binary and 
shock dissipation. The eccentricity of the disk in the final state 
of the disk varies approximately between 0. 1 and 0.16 depending 
on the position of the binary in its orbit as shown in the bottom 
panel of Fig. [3] In the outer parts of the disk beyond r ^ 5 AU, 
where the density is very low, the disk eccentricity is measured 
to be larger due to the proximity to the secondary star. 

The overall (mass averaged) eccentricity of the disk dis- 
plays oscillations on longer time scales of about 17 binary orbits 
around the mean value of 0.12 (left panel of Fig. |4|. The disk 
eccentricity Cdiskir) in that plot has been obtained by calculating 
the eccentricity of each disk element, as if in a two body mo- 
tion with the primary star, and then averaged azimuthally over 
the respective annulus. The superposed spikes are caused by the 
individual periapses of each orbit, during which the whole disk 
is perturbed by the strong spiral arms. At the same time the disk 
as a whole precesses, as is shown in the right panel of Fig. |4l 
This coherent slow precession is retrograde with a pattern speed 
much smaller than the orbital period of the disk material around 
the star. This interesting global disk phenomenon appears to be 
related to the induced precession of accretion disks in close bi- 
nary systems such as cataclysmic variable stars which occurs 
even in the case of circular binaries (|Lubow"19 91|). Very simila r 
results have recently been foun d bvlPaardekoo per et alj ( l2008l) . 
see also iKlev & Ne lson' (iQQ't). In test simulations using circu- 
lar and non-circular binaries we have also verified that eccentric 
disks occur in the circul a r case as well, and in recent simula- 
tions iPaardekooper et al] (|2008|) have also demonstrated (for y 
Cep parameter) that in the limit of circular orbits the mechanism 
that causes the instabili ty is indeed the 3:1 eccentric Lindblad 
resonance ( lLubowlll991 ). In a related study we have analysed 
in detail the influence of numerical parameter on the occurrence 
and properties of eccentric disks in binaries and we refer the 
reader to that paper for more details (Kley, Papaloizou & Ogilvie 
to be submitted). For the highly truncated disk in this study we 
expect the averaged eccentricity to be relatively robust against 
numerical variations. Test simulations of this initial equilibra- 
tion phase of the disk with different grid resolutions resulted 
in very similar values of the disk eccentricity. For the typical 
values of disk viscosity and temperature used by us we expect 
these values for the disk eccentricity, an expectation which is 
confirmed by Paardekooper et all ( l2008l) who find find very sim- 
ilar results. Concerning the evolution of embedded protoplanets 
in these disks it appears that the effect of the periodic disturbance 
of the binary may play a more important role than the eccentric- 
ity of the disk. 

The direction of the precession (prograde or retrograde with 
respect to the disk and binary motion) is influenced by pressure 
effects in the disk. For our chosen value of H/r - 0.05 we do in- 
deed expect a retrograde precession. Similar behavi our has also 
been demonstrated for disks with free eccentricity dPapaloizoul 
2005). 



3.2. The orbital elements of the binary 

In the previous section we have seen that the gravitational 
torques of the binary lead to truncation of the disk and re- 
arrangement of the material within. In turn, we expect a change 
in the orbital elements of the binary. 

To estimate theoretically the magnitude of the back-reaction 
a circumstellar disk has on the orbital elements of the binary 
we assume an idealised system consisting of a binary system 
and a ring-like mass distribution orbiting star 1 with mass m„„g. 
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Fig. 4. The evolution of the global (density weighted) averaged 
disk eccentricity (left) and the position angle of the disk's peri- 
apse (right). 

at a distance ((5-function) of rri„g. The energy Ebin and angular 
momentum Lbin of the binary is given by 

GMu. / 9 \i/2 

Etin = -i^, Lhi„^n[GMabi„{\-el„)) , (4) 

and the corresponding quantities of the ring are 

GMiMdisk _ /„ \i/2 

^ring — ^ 5 ^fing — faring I vj^*^l^r/ng I 5 

^^ring 

where M - M\ + Mi is the total mass of the two stars and 
- M1M2/M is the reduced mass. Now, suppose that the ring 
is shifted from its initial position r^.^^^ to a smaller radius r^.^^^ 
keeping all its mass. This radius change mimics the truncation 
of disk by the binary. Through this process the ring's energy and 
angular momentum are reduced from E". and L". to E^. and 

ring ring ring 

L^.^^^. By conservation of total energy and angular momentum 

E — E,'ifjg + Ejjifj E — Efif^g + Ebi^ . (6) 

We can calculate the corresponding change in the orbital ele- 
ments of the binary from Ef. and Lf. to E^. and if. . For the 

Dm bin bin bin 

binary parameter masses M\ - I.6M0, M2 - QAMq with ini- 
tial orbital elements aZ- - 18.5AU and - 0.36 we find for 

bin bin 

the shift of a ring with m„„g = 4 x lO^^M© and initial radius 
r" - 4.0 AU to a final radius of /. - 2.0 AU that the binary 

ring ring J 

elements change to o^.^ = 19.4 AU and e^.^ = 0.41. A quite sub- 
stantial change considering the smallness of the ring's mass in 
comparison to the stellar masses. But the closeness of the ring 
to the which is primary allows it to gain a substantial amount 
of binding energy. The calculation is approximate in the sense 
that the energy and angular momentum of the ring are calculated 
with respect to the non-inertial coordinate frame centred on the 
primary. 

We can now compare this estimate with the previous hydro- 
dynamical simulations, and plot in Fig. |5]the evolution of a/,,,, 
and Cbin for about the first 100 binary periods without a planet 
included. As demonstrated above, the binary expands as it gains 
energy from the truncated disk and increases its eccentricity. 
The increase in Cbhi does not lead to a decrease in the angular 
momentum however, since it increases its separation, see Eq. |4] 
Whenever the binary is near periastron the gravitational interac- 
tion with the disk is maximal which results in the strong periodic 
spikes in the binary elements. The change in the orbital elements 
of the binary is somewhat smaller than the estimated values be- 
cause /) the disk mass is smaller in the hydrodynamic calculation 
and //) disk mass and angular momentum are stripped off by the 
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Fig. 5. The evolution of the binary elements due to the interac- 
tion with the circumstellar disk around the primary star, without 
an embedded planet. One binary orbit refers to approximately 
57yrs. Left: ai,„(f); Right: ew„(f). 

secondary and are lost through the outer boundary of the com- 
putational domain. The loss through the (open) inner boundary 
of the disk is only marginal. 

4. Evolution with an embedded planet 

In the previous section we have analysed the change of the or- 
bital elements of the binary caused by the mass rearrangement in 
the disk. Now we turn to the embedded planet in the disk, which 
is inserted into the disk after an equilibration time of 100 binary 
orbits (nearly 6000 yrs). This rather time consuming procedure 
to generate the initial state is necessary to obtain realistic initial 
conditions for the growing protoplanet. At the time of insertion 
of the planet the remaining disk mass is rescaled (through So) 
to contain 3 Mj^p within the computational domain, while keep- 
ing its two-dimensional distribution unchanged. To follow the 
growth and evolution from planetary cores to massive planets 
the mass of the inserted protoplanet is 36MEarih orbiting the pri- 
mary star. This mass was chosen as it corresponds closely to the 
mass above which forming g as giant planets u ndergo the final 
phase of rapid gas accretion dPapaloizou & Ne lson 2005). The 
typical starting values of the semi-major axis and eccentricity 
are Qj, - 2.5AU and - 0.0. 

4.1. Exploratory models and numerical tests 

After inserting the protoplanet on a circular orbit at 2.5 AU we 
expect that its orbital elements will change due to the gravita- 
tional influence of the binary and the distorted disk. To differen- 
tiate the different contributions we decided to check the origin of 
the dynamical behaviour, through a variation of physical condi- 
tions. The standard model resembles the true physical situation 
where the planet feels the full influence of the binary and the 
disk which is perturbed by the binary. In the other set-ups we 
switch the various contributions on and off. Avoiding additional 
complications in this analysis, the protoplanet is not allowed to 
accrete here (face - 0) and its mass remains fixed at 36MEaiiii- 

The results are displayed in Fig. |6] where the semi-major 
axis and eccentricity of the planet are shown for four differ- 
ent set-ups. The standard model (h26, dark solid line) refers to 
the full model (including binary and disk) as just described, in 
model (h26a, dashed bright line) the mass of the secondary has 
been switched off to test its influence, and the planet evolves in 
the initially eccentric disk which becomes more circular during 
the evolution because of the absence of the secondary. In model 
(h26b, dark dotted line) in addition to having no secondary the 



Kley & Nelson: Planet Formation in y Cep 



7 




Fig. 6. Time evolution of the semi-major axis and eccentricity of 
the embedded planet for different physical set-ups, neglecting 
either the secondary (h26a,h26b) or the disk (h26c), see text for 
details. 



density in the disk has been azimuthally symmetrized keeping 
the radial distribution intact. Hence, this model suits as a ref- 
erence of what happens in the single star case. In the last model 
(h26c, fine dotted line) the secondary is present but the disk mass 
has been reduced, such that effectively only a 3-body problem 
is solved. The curves with the periodic bumps in Fig. prefer to 
the cases (h26,h26c) including the secondary while the smoother 
curves describe the situation where the secondary star has been 
excluded (h26a,h26b). The bumps occur with the binary period 
and indicate the 'kick' the planet experiences due to the inter- 
action with the spiral arms in the disk and the direct interaction 
with the binary companion when it is at periastron. The results 
show that the main contribution to the initial growth of planetary 
eccentricity Cp is in fact the eccentric disk. The eccentricity 
for models (h26 and h26a) rises initially with the same speed but 
due to the fact that, without the secondary star, the disk slowly 
circularises, the model (h26a) falls behind the full model (which 
reaches - 0.2 at f ^ 500) and ep begins to drop off after 
a time of about 400 yrs after insertion of the protoplanet. The 
maximum eccentricity reached for this reduced case (with no 
binary) is only ep = 0.1. In the symmetrized disk case with- 
out secondary (h26b) the planet behaves as expected for single 
star disks, it displays inward migration (the fastest for all cases) 
and the eccentricity remains zero. The case with only the binary 
does not show any change in the semi-major axis Op but rather 
a slow rise in ep. Additional test cases with similar set-up but 
different starting distances Op = 3.0 and 3.5 AU display similar 
behaviour For the largest initial distance 3.5 AU there are clear 



Fig. 7. Time evolution of the semi-major axis and eccentricity 
of the embedded planet for different physical set-ups. Here, the 
standard case (h26), a modified model where the dynamical in- 
fluence of the secondary on the planet has been switched off 
(h26al), and a model (h26tl) with a highly reduced planet mass 
(test particle with 5 x IQ^^^MEanh) are displayed. 



jumps in ep evolution visible at each individual periapse. Hence, 
the rise of planetary eccentricity above the disk eccentricity for 
the full model is a combination of disk and binary effects. 

To further study the consequences of the individual contri- 
butions we have displayed in Fig. Qin addition to the standard 
model a model (h26al) where just the dynamical (gravitational) 
effect of the secondary on the planet has been taken out. In con- 
trast to the previous model (h26a) the secondary still acts here on 
the disk. In addition we display a model (h26tl) with a highly re- 
duced planet mass of 5 x lO^^^MEanh- The migration rate of the 
first two models is very similar while reduced (nearly zero) for 
the latter due to the smaller mass. The eccentricity evolution is 
initially very similar for all 3 models, while at later times the 
test particle simulation displays the largest planetary eccentric- 
ity growth and the model without the secondary the smallest. We 
suspect that the faster growth of epia„eAt) for the test particle is a 
result of the vanishing migration and consequently larger influ- 
ence of the binary. In an additional plot (Fig.[8]l the evolution of 
the averaged disk eccentricity for these models is displayed. In 
the very small test particle case (model h26tl) the disk eccentric- 
ity is only influenced by the binary and the ejiskit) evolution is a 
direct continuation of the initialisation process (refer to Fig.|4]i. 
From Fig. [T] and [8] it is apparent that there exists a phaseshift 
between the eccentricity of the disk and that of the planet. As 
demonstrated above (cf. Fig. |4|i the variation of ejiskif) is cor- 
related with the precession rate of the disk, and is due largely 
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Fig. 8. Time evolution of the averaged disk eccentricity for three 
different physical set-ups. The standard case (h26), a modified 
model where the dynamical influence of the secondary on the 
planet has been switched off (h26al), and a model (h26tl) with 
a highly reduced planet mass (test particle with 5 x 10 '^M^art/i) 
are displayed. 



to the secular disk-binary interaction. The planetary eccentric- 
ity on the other hand is driven primarily by secular interaction 
with the disk (as shown above), and the phaseshift between disk 
and planet eccentricity is caused by the shift in disk and planet 
longitudes of pericentre induced by their relative precession (see 
Fig.[II]i. 

The two other models demonstrate that the disk eccentricity 
is clearly influenced by the presence of even a small mass planet 
of 36MEanh- In the case where the influence of the secondary 
has been switched off the planet's eccentricity is smaller as is 
the disk eccentricity. These results demonstrate that the com- 
bined evolution of planet and disk properties are rather intri- 
cately related in the case of binary stars. Hence, the evolution 
of the planet in this situation is more complicated in contrast to 
a pl anet embedded in an eccentric disk without a secondary star 
(see lPapaloizoull2002l) . and not directly comparable. 

To test the influence of numerical issues we perform for the 
reference case (h26) additional simulations with different set- 
ups. Among the purely numerical parameters altered are: the 
time step length as given through the Courant number /c, the 
force cutoff" radius r, as described in Sect. 12.41 and the grid res- 
olution. The reference case in this section has been calculated 
with fc = 0.5, r, = 1.0 and a resolution of 300 x 300. As 
displayed in Fig. |9] a change in the time-step and force cut- 
off to = 0.75, r, = 0.8 (dark dotted curve until t = 3300) 
has only a marginal influence. A change of resolution, here to 
NrX - 500 X 500 leads to a faster migration rate (30-40%) 
and slightly larger eccentricities. The qualitative behaviour (in- 
ward migration and saturation of eccentricity) is not changed, 
however. Given these encouraging results, we use for computa- 
tional reasons the reference resolution of 300 x 300 and fc - 
0.75, r, - 0.8 in our subsequent models. 

4.2. Models without accretion onto pianet 

Planetary cores form in the outer cooler regions of protoplane- 
tary disks beyond the so-called ice-line. In a binary star system 
the outer disk is affected most by the secondary, and to find pos- 
sible restrictions on the planet forming regions in the disk it is 
important to analyse the evolution of cores near the outer parts 
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Fig. 9. Time evolution of the semi-major axis and eccentricity 
of the embedded planet for different numerical set-ups. For the 
reference case (h26) we used for the Courant number and force 
cutoff /c - 0.5 and r, = 1.0, respectively. The second model 
(h28: dark dashed line, run only until t x 700yrs) used fc = 
0.75, r, = 0.8, which has also been used for the last model (h28f: 
light dashed) that has a grid resolution of 500 x 500. 



of the disk. To study the effect of initial position we start our 
embryos at 3 different locations in the disk 2.5, 3.0 and 3.5 AU 
always on a circular orbit, and choose again non-accreting cores. 
Because the initial characteristic growth time of the cores may be 
long, even in comparison to the orbital period of the binary these 
sets of runs constitute a test suite to estimate the orbital evolution 
of protoplanets in the disk. The results for the semi-major axis 
and eccentricity evolution of the planet are displayed in Fig. [10] 
where the only difference in the three cases is the release distance 
(2.5, 3.0 and 3.5 AU) of the planet. From all three locations the 
planets migrate inwards at approximately the same rate with the 
tendency for a slowdown for the two outer cases. This inward 
migration of the outer planets which are mainly orbiting out- 
side the disk is due the the planet disk interaction near planetary 
periastron where the planet moves faster through the disk and in- 
duces a negative torque and power (energy loss) (ICressweU et al.l 
|2007|) . However, the different initial starting radii lead to a very 
different eccentricity evolution. While only the innermost planet 
(starting at 2.5 AU) shows a finite eccentricity evolution the two 
outer cases display a very strong increase in their eccentricity be- 
yond gp - 0.5 after about 55 binary orbits. Clearly the strongly 
disturbed disk in the outer regions at around 4 AU significantly 
perturbs the orbits of the protoplanet and does not allow for small 
equilibrium eccentricity. While the initial faster rise of e,, ap- 
pears to be caused primarily by the disk, the long-term evolution 
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Fig. 10. Time evolution of the semi-major axis and eccentric- 
ity of a non-accreting planet embedded planet for three different 
starting locations. 



to very high eccentricities is most likely due to the combined ac- 
tion of the binary and the eccentric disk. A test simulation (for an 
initial semi-major axis Oj, -3.0 AU) where the gravitational ac- 
tion on the planet due to the secondary has been switched off did 
not display any large eccentricity growth but a steady decline. On 
the other hand, starting the planets with high eccentricity (from 
Gp -3.0 AU) in a pure 3-body problem does also not lead nec- 
essarily to a growth of eccentricity but rather to an oscillatory 
behaviour. So we may conclude that in the present situation it is 
the combined action of an eccentric, periodically perturbed disk 
and the influence of an eccentric binary star. 

As the planets move on non-circular orbits in an eccentric 
disk and binary a temporal change of the apsidal line may be ex- 
pected. In Fig. [TT]the evolution of the argument of pericentre of 
the planets and the disk are plotted versus time. The innermost 
planet has on average a periastron angle of about 200 deg while 
the outer planets have a larger angle. For the binary the angle of 
periastron lies fixed at degrees, and the disk is slowly precess- 
ing retrograde as pointed out above and indicated in Fig. [TTIbv 
the tightly dashed line. The influence of the disk on the planets is 
visible by the slight oscillations of the periastron angle about the 
mean with the same period as the oscillations in the eccentricity. 
The innermost planet has approximately a phase shift of 180 deg 
with respect to the binary and is nearly in an anti-symmetric state 
while the other planets are lagging behind this configuration. 




500 1000 



1500 2000 
Time [Yrs] 



2500 3000 3500 



Fig. 11. Time evolution of the argument of pericentre (perias- 
tron) for an embedded planet starting at 3 different locations and 
the protoplanetary disk. 



4.3. Models with mass accretion onto pianet 

To make the scenario more realistic we add now the option of 
mass accretion onto the protoplanet. To do so, we assume that 
initially the planets grew during a long time to their present mass 
and enter now the more rapid gas accretion phase. Numerically 
we continue the previous simulations for the three initial radii 
at a time t = 2500 yrs. The evolution of the orbital elements of 
the planets are displayed in Fig. [12] where the time offset has 
been reset after the turn-on of mass accretion. The thinner lines 
(until t a: 3000) are the models with accretion switched on at a 
low rate using an efficiency of face = 3 x 10"^. The symbols re- 
fer to the final phase of the models without accretion. While the 
eccentricities of the outer planets continue to rise to very high 
values despite a slow inward migration, the inner planet's orbit 
(starting at Op - 2.3AU) begins to circularise while the inward 
migration speed increases. The large eccentricities of the outer 
planets can probably only be maintained in a stable way because 
of some eccentricity damping action of the disk during their peri- 
astron passage. The planetary mass evolution for the three cases 
is displayed in Fig. [13] The innermost planet has reached about 
half a Jupiter mass after 3400 yrs while the other two grow at 
a smaller rate and have reached about 0.3 My„p after 3000 yrs. 
Note that these labels refer to the initial positions of the planets 
that were first evolved without mass accretion and then contin- 
ued with face ^ as described above. From the simulations it 
seems to be clear that it may be difficult to grow a planet at too 
large a distance from the primary star. For our chosen binary pa- 
rameter the borderline lies around 2.7 AU, which has been ver- 
ified by additional simulations with different starting distances 
that are not shown here. This particular transition radius of about 
2.7 AU depends on the chosen parameters for the disk viscosity 
and temperature. Larger values can possibly extend this transi- 
tion radius slightly while smaller values will reduce it. Varying 
the disk mass will change the viscosity and temperature of the 
disk, and we expect this to influence the transition radius, too. It 
is beyond the scope of the present to paper to perform detailed 
parameter studies on these issues. From the results displayed in 
Fig. [13] it appears that mass accretion onto the planet does not 
alter the orbital evolution of the planet significantly. 

To further explore parameter space we used the starting con- 
dition Gp = 2.5 and performed one additional run with a three 
times higher accretion rate onto the planet (face - 10 denoted 
'medium' in the figures). The results of this new model with 
the previous no accretion and low accretion models are shown 
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Fig. 12. Time evolution of the semi-major axis and eccentricity 
for 3 accreting protoplanets starting at 3 different locations in the 
disk. The symbols refer to the non-accreting models of Fig. \W\ 
time-shifted by 2500 yrs, while the lines refer to restarted mod- 
els with accretion switched on. 
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Fig. 13. Mass of the planet vs. time for accreting models for 3 
different initial locations of the planet. 



in Fig. [14] The results for the orbital elements of the planet 
in Fig. [14] show that the models with no and low accretion be- 
have similarly. The planet continues to migrate inward while de- 
creasing its eccentricity. The model with the higher accretion 
rate migrates faster at first, and slows down near the end, while 
the eccentricity settles to about Cp = 0.05. This reduction in the 
migration rate is a consequence of the increasing mass of the 
planet. For this higher accretion rate the planet's mass increases 
to nearly 1.7 Mjup after 3300 years (see Fig. [15). Since the disk 



Fig. 14. Time evolution of the semi-major axis and eccentricity 
of 3 accreting protoplanets starting at the same initial condition 
with different mass accretion rates. 
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Fig. 15. Mass of the planet vs. time for accreting models starting 
at Op = 2.5AU for 3 different accretion rates. 



mass is reduced by the same amount the driving agent of the 
planet is lost and the speed of migration reduced. This model 
with the 'medium' accretion rate settles to a final state which is 
quite close to the observed value for the y Cep system. 

A massive embedded planet will open a gap in standard cir- 
cular disks, and it is interesting to analyse this effect within the 
present context. In Fig.[T6]we display the two-dimensional den- 
sity distribution l.{r,(p) in the disk at a time 3125 yrs for the 
medium accretion model. At this time the planet has reached 
approximately a mass of 1.5 Mjup. From the plot it seems that 
the disk inside the planetary orbit is apparently more circu- 



Kley & Nelson: Planet Formation in y Cep 



11 




-6 -4 -2 2 4 



X 

Fig. 16. Grayscale plot of the two dimensional density distribu- 
tion of the medium accretion model (h48) at time 3125 yrs. The 
shading is scaled oc E'^"* between 4.8 x 10"'' (white) and 2400 
g/cm^ (black). The location of the planet is marked by the small 
circle. 
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Fig. 17. Azimuthally averaged density (in g/cm^) and eccen- 
tricity of the disk for the medium accretion model at a time of 
3125 yrs. 

lar than outside. This is confirmed by the corresponding one- 
dimensional radial distribution of the azimuthally averaged den- 
sity and eccentricity of the disk at the same time. The gap is 
somewhat weaker than in circular disks primarily due to the pe- 
riodic disturbance of the secondary that tends to sweep material 
into the cleared region around the planet. Due to the shallower 
gap the planet is able to continue mass accretion from its sur- 
roundings more easily compared to a planet on a circular orbit in 
a single star system. The inner disk clearly has a smaller eccen- 
tricity than the outer parts (compare this to Fig. [3]). The presence 
of the planet represents a barrier for the (spiral) wave induced 
by the binary which consequently cannot propagate into the in- 
ner parts of the disk. 

5. Conclusions and Summary 

In the present work we have concentrated on the planetary for- 
mation in the system y Cep which places, due to its binary pa- 
rameter of Ohin = 20AU and eft,„ - 0.41, severe constraints on 
the formation process. We investigated an intermediate phase in 



the planet formation process within the so-called core formation 
scenario, i.e. we did not attempt to model the formation of plan- 
etesimals in this system but rather concentrated on a later phase 
when protoplanetary cores have already formed and begin their 
rapid gas accretion phase. 

Before embedding the protoplanets in the disk, we have first 
brought the system into equilibrium by performing simulations 
with no embedded protoplanets. This initialisation phase takes 
about 100 binary orbits after which the disk structure has set- 
tled to a new quasi-equilibrium configuration with a truncated 
disk. Interestingly, the accretion disk around the primary reaches 
an average eccentricity of about ejisk = 0.12 and shows a co- 
herent retrograde p r ecessi on, in agreement with the findings of 
iPaardekooper et all (|2008|) . A detailed parameter study of this 
interesting behaviour on its own is beyond the present analysis 
and will be presented elsewhere. During one binary orbit the ec- 
centric secondary disturbs the disk periodically and induces sig- 
nificant non-axisymmetric perturbations in the disk which decay 
due to viscosity during binary apoapse. 

As suspected these perturbations of the disk, in particular the 
periodic creation of strong tidally induced spiral density waves 
and the creation of an eccentric disk, lead to non-negligible effect 
on the planetary orbital elements. While embryos placed in the 
disk at different initial distances from the primary star continue 
to migrate inward at approximately the same rate, the eccentric- 
ity evolution is markedly different for the individual cases. If the 
initial distance is beyond about a^2J AU the eccentricity of the 
embryo continues to rise to very high values and, apparently only 
due to the damping action of the disk, the orbits remains bound. 
The excitation mechanism of the eccentricity is the combined ac- 
tion of the binary and the perturbed disk (i.e. there is a long term 
secular interaction due to the disk being eccentric). Test sim- 
ulations where individual components of the system have been 
switched off and on have made it clear that it is indeed the combi- 
nation of effects (eccentric binary and the eccentric, periodically 
perturbed disk) that leads to the high planetary eccentricities in 
the system. Once epia„et has grown to sufficiently large values 
the planetary orbit may reach the stability limit as given for ex- 
ample by the pur e 3-body simulations for y Cep (see Fig. 2 in 
iTurrini et al.ll20(j4i) . 

A low mass, non accreting planet embedded in an the ec- 
centric disk experienced substantial growth in eccentricity (see 
Fig. [Tol l. This has clear implications for the accretion of plan- 
etesimals because thei r velocity dispersion may become very 
large due to this effect. iThebault et akl (|2004|) examined the evo- 
lution of planetesimal orbits under the influence of the binary 
companion and aerodynamical gas drag. They concluded that 
accretion of planetesimals would occur in the shear dominated 
regime because orbital alignment was maintained due to the gas 
drag. This work, however, did not include the effects of an eccen- 
tric disk, and so it remains unclear whether planetesimal orbits 
will remain ahgned. R ecen tly, this issue has been add ressed by 
iKlev & Nelsonr(l2007l) and IPaardekooper et all (l2008l) who find 
that the inclusion of the dynamical evolution of the disk may 
render the planetesimal formation more difficult, but this topic 
certainly deserves more consideration. 

The most up to date observational data suggest the follow- 
ing parameters for the planet in the y Cep system: Qj, ^ 2.044, 
ep - 0.115 and nip sin i ^ 1.60 Mjupiwr- If this planet formed 
according to the core instability model, then an important issue 
is the survival of the planetary core before gas accretion occurs. 
Fig. [To] has shown that the non-accreting, low mass planet un- 
dergoes quite rapid inward migration. The migration, however, 
is modulated by the eccentricity of the planet, such that at high 
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eccentricity phases the migration rate decreases. It is possible 
that longer run times will show an essential stalling of this mi- 
gration if the planet eccentricity grows beyond its final value of 
Bp ^ 0.3. The inclusion of a better approximation to the physical 
state of the disk, eg. radiative cooling, may help to slow down 
this initial phase, as new results indicate that in case of radiative 
disks migration is slowed down or even reversed for protoplane- 
tary c ores (Paardekooper & Mellema 2008; Baruteau & Ma ssed 
120081) . 

Once gas accretion is switched on, it is clear that a disk 
mass of about 3 Jupiter masses, where the outer disk radius is 
tidally truncated at r ^ 5 AU, will be sufficient to grow a planet 
that is close to the minimum observed mass of nip sin / ^ 2.044 
Mjupiier- It is also clear that we can construct a model in which 
a low mass planet growing from an initially circular orbit can 
achieve a final mass of nip ^ 2 M jupiier, and have a final eccen- 
tricity of Cp ^0.1 as required. 

A final comment relates to the final mass of the planet. Our 
simulations suggest that a disk mass of about 3 Jupiter masses 
will be enough to form a gas giant of the required minimum 
mass. A future test of the mode by which the planet in y Cep 
formed (gravitational instability versus core accretion) will be 
determination of its actual mass. We suspect that a disk that is 
massive enough to form a planet through gravitational instability 
will lead to a planet whose final mass is substantially larger than 
the minimum value observed. 
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